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ABSTRACT 


A computer program for determining the theoretical performance of a three- 
dimensional inlet is presented. An analysis for determining the capture area, 
ram force, spillage force, and surface pressure force and their moment is pre- 
sented along with the necessary computer program. A sample calculation is also 
included. 
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INTRODUCTION 


In the design of an integrated hypersonic vehicle, the various forces 
generated by the airbrenthing propulsion system must be determined and located 
accurately with respect to the center of gravity to determine the interaction 
between the propulsive unit and the vehicle. 

The supersonic combustion ramjet engine (scramjet) has been suggested as a 
possible airbreathing propulsion system on a hypersonic vehicle such as a launch 
vehicle.^ The scramjet engine has three basic components: the inlet, combustion, 
and nozzle. For hypersonic flight these components tend to be of sizable 
dimensions with respect to the vehicle size. Therefore, the forces generated by 
the scramjet engine especially those of the inlet and nozzle are also large com- 
pared with other aerodynamic forces. Thus an accurate evaluation and location of 
these forces is essential to determine the trim requirements of the vehicle. 

In this report, the pressure and momentum integrals necessary to establish 

the external forces generated by the inlet are evaluated in conjunction with a 

2 

three-dimensional computer program which determines the external flow field 

3 

using the theory presented in . Although the flow field analysis is inviscid in 
nature, the body geometry can be corrected by adding a suitable boundary layer 
displacement thickness to the actual body. 


t 


ANALYSIS 


A typical scramjet engine installation is shown in Fig. 1. The propulsion 
system is surrounded by a suitable control volume to evaluate the integrals 
necessary in the determination of the propulsive forces. The control volume 
shotm in Fig. 1 assumes the captured streamline in the vehicle flow field at 
station "0". That is, the streamline crosses the bow shock upstream of this sta- 
tion. This is typical of low hypersonic speeds. At hypersonic speeds, the captured 
streamline may have crossed the bow shock at a station downstream of station "0" 
as the bow shock is very shallow. The control volume changes not only with 
flight Mach number, but also with vehicle attitude (i.e. angle of attack, yaw, 
etc.). Also, here we assume the inlet is started and there is no flow separation 

present on the inlet ramp. The particular inlet performance evaluated here are 

\ 

the local capture area, ram force, spillage and body surface, forces and moments. 

The forces are evaluated in terms of components in the axial, normal, and side 
wise directions of a body fixed coordinated system. 

a) Local Capture Area 

To determine the inlet's capture area or mass flow, it is necessary to de- 
termine the trace at station "0" of the stream surface that impinges on the inlets 
cowl and sidewalls leading edges. Since the program marches forward in the axial 
direction, it is not known a priori through which initial point the captured 
streamline will pass. Therefore, a number of streamlines must be traced from the 
initial station. The streamlines intersecting the leading edges are then selected 
or interpolated at the end of the calculation. The program has been set up to trace 
streamlines passing through the even numbered index points in the radial direction 
at the initial plane (i.e. J = 2, 4,..., K=l, 2, 3,... See Fig. 2). 


The trace of a generic streamline proceeds as follows. The equations of a 


streamline in cylindrical coordinates are 


dx dv 

s _ s 

U V 



where u, v, w, are the axial, radial and outflow components of velocity in 

cylindrical coordinates, and the subscript s denotes streamline. To trace a 

streamline which at the initial station (x = x^) passes through the grid point 

(r , 6 ) indexed by J,K proceed from the given initial data plane by estimating 
o o 

the position of the streamline at the next calculation step i.e. 
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where the superscript (1) denotes first iteration and the subscripts i, j, k- 
denote data at the present station x. and r g 

X Sf Si 

The three-dimensional program then proceeds to calculate new data at the 

next axial station x. = x. + £x. Given this data and the first estimate of the 

ItI x 

streamline position at the new station the end slope is determined from inter- 
polated data as explained in Appendix A. The average slope between station x^ 

and x. ' is used to determine the location and slopes of the streamline at the 
l+l 


new station. No additional iterations are taken as the calculated step size is 

sufficiently small that the flow properties do not change appreciably. The 

-4 

error has been found to be on the order of 10 

The computer program subroutine outlining this procedure is presented in 
Appendix B. A vertical plane of symmetry is assumed in the program at tho piano 
indexed by K = 1 . In this plane the program calculates only the strearraise slope 
dr/dx and sets q = 0. 

Additional streamlines are traced at a station downstream of the initial 
station to insure that the cowl leading edge is bracketed by the selected stream- 
lines. The additional streamlines originate at odd numbered index points in the 
radial direction (i.e. J= 1 , 3 , 4 , 7 ,... for K=l, 2 , 3 , . . . ) . The station ac which the 
additional streamlines are initiated is determined from the mean cowl lip radius 
and the conical shock angle by using an approximate formula for the streamlines 

through the conical flow. In the approximation the streamline is a hyperbola 
2 2 2 

(Ref. 4 ): r = (tan 0 c ) x + C where C is a constant. Selecting the constant 

so that the streamline impinges on the cowl lip i.e. @ x = x^ r = gives 

2 2 2 2 2 2 
r = (R - x tan 0 ) + x tan 0 . 
c c c c 

r = x tan 0 at: 
w 

x = x 

c 

where 0 is the conical shock angle 
w 

6 c is the cone angle 

and g is the ray angle passing through the cowl lip = Tan \r c /x c ) 


The streamline crosses the conical shock 


2 2 
tan 0 . - tan q 

cowl °c 


2 

tan 0 


w 


- tan 9 


The effect of spillage is to bend the streamlines further from the conical 
flow above. Therefore, this choice of initial station to trace the streamlines 
will insure that some of the selected streamlines will be spilled arcund the cowl. 
However, interpolation and/or extrapolation of these streamlines may still be 
necessary upstream o£ this station for odd-indexed captured streamlines that lie 
near the shock wave. 

Given the above set of streamlines, the boundaries of the local captured area 
are determined by interpolating between the traced streamlines from the cowl lead- 
ing edge to the initial station as described in Appendix A. 

b) Captured Local Mass Flow Rate 

Given the boundary of the local capture area at station "0", the local captured 
flow rate is determined by 


4 * JJ P" “o 


A' 

o 

where A^ is the local capture area at station "0" as shown in Fig. 1. 

However, since the limits of integration are not known at the start of the 
program and to keep' track of the accuracy of the streamlines tracing procedure, 
the incremental mass flow between mesh grids is calculated at every printout 
station,- i.e. 

Am = 0 u fck 


or 
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(P/P M) 


^ - J W+& 
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A . 
ref 


(!) 


P A _ 
oo ref 
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where *A ~ r A 6Ar = % (r j>k + r.^ + r j+1#fcH + ^[(r^^- r j>fc ) 


(r j+l, k+l“ r j,k+l^ 48 


The quantity in brackets is evaluated at the four corner points or each grid and 
the average 1 is printed out at every print -out station. The total captured mass 
flow is the summation of all grids in the captured streamtube, i.e. 
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i k 
J capt, capt 

j >k 

This quantity can be nondimensionalized with respect to the maximum local captured 
flow or can be converted into equivalent freestream tube capture area 
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O CO 03 


r 


Bi 


-M Y 
1 00 ' 


L 


J/(T/T o ) 


] 


O oo 


g. p 
a co 
o 


and nondimensionalized with respect to the maximum freestream capture area or cowl 
frontal area, A . 
c) Ram Force 

The ram force in general has three components: axial, normal, and side. For 
zero yaw, the third component is zero due to symmetry. The individual components 
are calculated from the incremental ram force in the grid and then summed as in 
the mass flow calculation. 

i) axial component 

the incremental axial component of the ram force is given by 


h * 1 


ram 


axial 


= A*u + (p~P ) A a 


or &F 


ram 


axial 
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AA 


where AA is given above. 


The computer program calculates this quantity at the four corners of the grid and 
prints out the average. 

The total axial component of the ran force is the double sum of all grids in 
the captured streamtube. The axial component ram force coefficient is then given 
by 
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ram 

ax 
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ref 


ii) norma 1 - c omp onent 

The incremental component of the ran force normal to a plane perpendi- 
cular to the vehicles plane of symmetry is given by 

* F ram , = **“<5 cos 9 ' Z Sin 6> 
normal 

AA 


iii) side force component 

The incremental side force component is given by 

aF = Amu (— sin e + — cos e ) 

*- ram u \u u °/ 
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The calculation procedure to obtain the normal and side ran force components 
coefficient is identical to the above. The tota'l ran force side component 
vanishes for zero yaw due to symmetry requirements, 
iv) ram force direction 

The angle which the ram force makes with the vehicle's axis is given by 


e -Tan- 1 /^ L N, 

\ ram . / ram . . / 
normal / axial 

and with the y-coordinate is given by 



v) moment produced by the ram force 

The incremental pitching moment produced by the ram force about any 
point such as the c.g. is given by 



■where the subscript q denotes a pitching moment about the y-axis. 
The incremental rolling moment is given by: 

K ) - (^*) ac f - (^*) 4C. 


normal 


ram . , 
side 


where x,y,z are the coordinates of the centroid of the incremental area. 


The total pitching, rolling and yawing momentsare calculated by sunning over 


all grids inside the captured streamtube of the inlet, that is 
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^ j,k 
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For a vehicle with a vertical plane of symmetry and at zero yaw only' a pitching 
moment is produced. The computer program calculates only the pitching moment, 
about the origin (cone tip). The distance to the center of area is used in the 
calculation. The calculation is made at every print out station 

vi) Point of Application 

The point of application of the ram force is where the line of action of 
.the ram force pierces the y-z plane at station "0". This point is found by 
locating the point where the ram force generates the same pitching, rolling and 
yawing moments , that is : 
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Only two of the above equations are linearly independent. The third expresses 
the orthogonality condition between the moment vector and the force vector, i.e. 
M®F = 0. Therefore only the line of action of the ram force can be found from 
the above equations, i.e. _ _ _ 
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for a vehicle with a vertical plane of symmetry and at zero yaw, the second 
equation is indeterminate. Further, selecting x = x q gives the piercing point 

fsj 

z. 

d) Surface Pressure Integral 

The pressure integral on the inlet ramp has three components; axial, normal, 
and side components. Since the limits of integration are not knoxm initially, 
i.e. the limiting streamline impinging on the sidewall leading edge intersection 
point with the surface (point "A" in Fig. 1), the incremental surface pressure 
integrals are calculated at each calculation step, 
i) axial component 

The incremental axial component of the inlet ramp is given by 


^ ^ax % ^axial^ ^ref 


AS ax « r A9 &r 


where C is the average pressure coefficient (C =(P-P )/M %v P ) of the C 's at the 

P P OS' OD 1 OQ P 

four corners of the grid, and AS ax ^ a ^ is the incremental body surface area in a 
typical grid (Fig. 2) AS. The total axial force is formed by the summation 
over all grids from the station "0" to the inlet's "mouth", 
ii) normal component 

The normal component is calculated from 


AC„ = C AS £ 

• u N p N rer 


AS N % (Ax) (Ay) 


where is the incremental body surface area an a typical grid projected onto 
a horizonaal plane normal to the plane of symmetry- of the vehicle 

iii) side component 

The side component is calculated from 
iC s * C p ‘ S s /S ref ^ S s * to) to) 

where AS & is the incremental body surface projected onto the vehicle's plane of 
plane of symmetry. 

iv) direction of the body surface 

The body surface force direction is given by 

? - to ' 1 <c K /c ax ) 

and 

C - W 1 (C s /C ax ) 


v) moment produced by body surface force 

The incremental pitching moment produced by the body surface force about 
any point is given by 
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The incremental rolling moment 
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The incremental yawing moment 
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Total pitching, rolling, and yawing moments are obtained by sunning over all the 
grids inside the welted surface from station "0" to the "mouth" of the inlet. 

For a vehicle with a vertical plane of symmetry and at zero yaw, the 
rolling and yawing moments vanish. The computer program calculates only the 
pitching moment at every station vhere printout ggguis. 

vi) point of application of body surface force 

The point of application of the body pressure force is where the line o 
action of the body surface force pierces a horizontal plane (z=0). This point 
is given by locating the point where the body force produces the same pitching, 
rolling and yawing moments , that is : 

Si ‘ 2 ' C AX ' J °N ( + 
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As before, only the line of action can be found 
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For a vehicle with a vertical plane of symmetry and at zero yaw the second 
equation is indeterminate. Further, selecting z=0 gives the piercing point x. 


e) Spillage Force 


The spillage force components can be calculated in analogous manner as the 
body surface force except that the spillage surface and the flow properties on it 
are used. However, this surface and its properties are not knotm a priori. 
Therefore it is more convenient to determine this force by imposing the equilibrium 
of a control volume consisting of the ram force at station "0", the body surface 
force and the ram force at the mouth of the inlet. The resultant of these is the 
spillage force which includes the side wall and cowl spillage. Equilibrium of 
moments yields the point of application and the moment of the spillage force. 

That is, the spillage force is given by a vectorial sum of the forces acting on 
the above volume: 

-F = F + F + F 

spill ram„ 0 „ body ra ™ cdwl 

from which the three scalar components can be obtained. For a vehicle at zero 
yaw and with a vertical plane of symmetry only the axial and normal components 
are nonvanishing and are equal to: 
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The pitching moments produced by the spillage force is givenby 
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The line of action of the spillage force is given by 
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The angle which the spillage force makes with the vehicle's axis is given by 


8 =■ Ian" 1 /C, 


'spUy 


' spill 
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and with the y-coordinate is given by 


T1 " Tan ~ 1 f c v 
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SAMPLE CALCULATION 


A sample calculation of the above forces was done for the inlet geometry 
described in Ref. 5. The free stream Mach number was equal to 6.0 and the 
vehicle was at zero angle of attack. 

The streamlines in the vehicle's plane of symmetry and the projection of the 
streamlines onto a plane 6“ to the vehicle axis are shown in Figs. 3 and, 4 
respectively. 

The captured streamtube at station "0" is shown in Fig. 5. The projections 
of the limiting streamlines are also drawn. The captured mass flow, ram force, 
body force, and spillage surface, directions and magnitude as well as moments are 
listed in Table 1 for this condition. 


RESULTS 


A computer program for evaluating the captured mass flow, ram, surface, and 
spillage forces as well as their moments, and points of application for a three- 
dimensional inlet design have been presented along with a sample calculation. 

The usefulness of this tool in determining the inlet performance and forces 
and moments produced by the inlet's ramp surface is appreciated in a trim 
analysis of a hypersonic integrated vehicle. 
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APPENDIX A 


1) Interpolation Scheme for Properties at a Point Inside a Grid . After 
locating a point inside a mesh whose corner points are indexed by (j,k), (j+1, k), 
(j, k+1), (j+1, k+1), (e 1 in Fig. 2c) the properties at the interior point are 
evaluated through a linear interpolation 
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Except at the plane of symmetry where the interpolation is carried out with 
respect to r only, i.e. 
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while for the surface streamlines (j = 1, k), the interpolation is carried out , 
with respect to 0 only, i.e. 
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2) Interpolation Scheme for Captured Streamlines Given the location of the 
cowl leading edge piercing point (e* in Fig. 2c) with a planh x = constant, four 
adjacent streamlines (aa 1 , bb', cc', dd') which enclose the point e' are used in the 
following interpolation scheme -to determine the location e at the initial station 
of the streamline passing through e 1 . 


For any general poxnt, the location of the point e xs given by 
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, In the present tracing scheme 
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For streamlines in the plane of symmetry a simple interpolation is used. 
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APPENDIX B 


PROGRAM LISTING 


* 


•k 

* 

* 

* 

** 


* 


* 


* 


PROGRAM GRUM( input, output , punch, TAPE5 = INPUT,TaP_E 6=0UTPUT, 
lTAPE7 = PU.NiCn»TAPE8) 

DIMENSION Ul(21,20) tV 1(21,20) » Wl ( 21, 20 ) , 

1R2(20) »RX2(20)»RT2(20>%DR2(20) , DRT2(20) ,RS2(20) ,RSX2(2D) , RST2<20) , 
2F21(21»20»b) ,F22< 21,20,5) ,F23(2i« 20,5) , 

2U2(21«2Q) , V 2(21*20) , W 2(21, 20), EL 13(21, 20) , 

3ELP3 ( 21,20) ,Uf 2 (21,20 > , VT2 ( 21 « 20 ) , 'JT2 ( 2l , 20 ) , £L2 { 21 , 20 , 2 ) , , 
4R3X(2),RA(2),fA(2,E) , EL A { 2 , 2 ) , GA { 2 ) , G3A ( 2 ) , uEPEP ( 2 ) ,CA?A{2) , 
•5S5AR(2) ,ELfcJAR(2) , EM (2 1,20) «CP(2l«20) 

DIMENSION XpROBE { / ) ,°3AR(21,20) 

COMMON /FIK/X3(20) ,RB( 20,20) »R3P(20,2Q) , R3PP ( 20 , 20 ) , RB? 33 ( 2,0 , 20 ) , 
luu( 20) «VU(20> »1.'U( 20) ,R1 (20) ,RX1 (20 ) ,PT1 (20 ) ,DRl(2Q ) ,QRT1 (20 ) , 

2RS1 ( 20 ) ,RSX1( 20) ,RST1 (20) *Fll{ 21,20,6) ,F12{21, 20,6) , Flo (21,20, G) , 
2F14(21,20,M ,EL1( 21,20,2) , T ( 20 ) , P ( 21 ) , 

3 T 0 LEP»CVNST,crST, CN,v , CNN* , xEftIO , 0 T , DTD , XI , VI WF2 , 52 , G3 , G4 , 35 , CCPl , 

4 C C P 2 , M I , F J , H K , M A X M , M A X N , M 0 D E , I P T , K P T , M 1 5 M 1 , M J M 1 
COMMON /F2/S(30)/r‘l2/FPPP(5Q)/Fi3/0S/Fl23/F(50) , FP ( 50 ) « F 3 P ( 50 ) 
1/C1/G1, AO/riESH/NOPT, I 3 UNCH 

E3UI VALENCE ( F1K1 ) ,U1(D ) , (F1K421 ) , Ml ( 1 ) ) , ( Fll ( B4l ) . •„ 1 { 1 ) ) 
E0UIVALENC£(F21(1) ,u2(l) ) , (F22(l ) , Ut2 ( 1 ) ) , ( F21 ( 421 ),V2(1)),{ 

1 F 2 2 ( 4 2 1 ) , VT2(1) ) , (F21(641) , W2(1) ) , (F22{841) ,wT2<l) ) , ( 

2F2l( 1261 ) ,£li 3(D > , ( F21 ( 1&31 > , EL 23 ( i ) > 

COMMON/STLIN/G, X2,RR, GRQX ,D 70X , RrS A V , THSAV. EM,CP,P3AR 
. DIMENSION URDX(21»20) , DTQX < 21 , 20 } , RPSAV ( 21 , 20 ) , THS AV ( 21 , 20 ) , 

1RR ( 21 , 20 ) 

ICYCLE=0 

MCYCLt=0 

MCOWvJ=0 

C 

C MC0WW=0 TRACE only even numbered streamlines 
c M C 0 W W = 1 TRACE EVERY STREAMLINE 

• REAQ(5.1) lMINF , (XP.R03E( I) , 1=1,6) , TCOWL 

1 FORMATt 6F10.5) 

XC=11, 

XPRINT=XPR03E(1) 

IPROBE = 1 
ICO = 0 
I JL=0 
I W = 6 

call FIRST 

XCOWW=XC*SQRT( (TAM (TCOWL) **2-TAN ( RBP ( 1 , 1 ) ) **2 ) / ( TAN ( RSX1 ( 1 ) ) **2- 
1 TAN(R3P(1, 1) >**2) ) 

• WRITE (6,2) XCUWW 

2 FORMAT!///, 5X,* 00D NUMBERED STREAMLINES WILL BE STARTED AT X= *, 

1 E15 , 5/ ) 

G= G4/(G4-1.) 

So = 2 « / ( G-l , ) 

CONST = .5*G*EMINF**2. 

1 = 0 

GO TO 2201 " 

C 

C SUBSEQUENT INITIAL OATA SURFACE (2100) 

2101 Xl=X2 

DO 2104 K = 1 , MR 

— • Rl(K)=R2(K) ‘ ~ ' 

RX1(K)=RX2(K) 

RT1(K)=RT2(i<) 



0 

* 

c i 

O l 

03 6006 J = 1,."U - ‘ 

R?=RP+0R2(K) 

* RR ( J i K ) =RP 

V?2 = U2( J,‘<)#*2+\/2( J,K)**2+W2( J«K>**2 

A=A0-Gi*VP2 

ZVZ-\l?2/k 

E l U J,K)=S3KT( EM2) 

RiTA=S0RT(L,M2-l. ) 

TANALF=ABS( LJ2( J«K) ) /SORT < U2 C J. K ) *#2- V2 < J , K > ) 

C L A S 3 - ( B C T A T A Ns A L “ + 1 . ) / ( b £ f A - T A f>J A L F > 

TcST=OX/DT-CTST^{ RP*U2< J,K)-DX*v2< J,K) ) / ( ELAM3*SQRT ( U2 C J » K ) **2 + 
1V2( J«X)4.*2> ) 

IF ( TEST > 6005, 6005, 9.002 
&005 CPC J,K) = (PIRAr-i=A<-«-GA-CCPi>-'-CCP2 
Y = 'DR2(K)*FL0AT( J-l) 
t v 112=EM( J»*)**2. 

P 3 A R ( «J » K ) = 1. +C 3 i J,<)*COfJST 

PT23AR =PTKATfc(E‘-U2/G2/ (EN12 + G&) ) **G4* ( G2* C 2 . *G«*#EM12-1 . ) ) *s(l./ 
l(l.-G) ) 

GO TO (6008,6006) N?R T 

6000 »JRITE (IW* 20021} J , RP « U2 ( J , K ) « V2 ( J« K ) , W2 ( J, K ) » EL2 C J, K . 1 ) . 

1EL2( J,K,2) .EMC J , K ) , CP < J , K ) , Y . P3A R ( J , K ) ,PT23AR 
FORMAT (I5,11F10,6) 

CONTINUE 
CONTINUE 
ICYCl£=ICYClE-H 
I F ( X2 .LT. XCOwi-J) GO TO 12 
MC0*IW = 1 

v iCYCLE = MCYCL£*l 
CONTINUE 
GO TO (8.6) NPRT 
/J R I T E ( 6 , 9 ) X2 

FORMAT ( /// » 5X« ♦STREAMLINE TR ACE* , 5X , *X2 =*,E12.5) 

CONTINUE 

CALL STREA'HAPRT, V 1 J, M<« I C YCLE. MC YCLE . MCOWW . OX , XCOWH. U2 , VZ . W2 . T ) 
IFCX2.GE.5.0) N0PT = 4- 
& 

d 
C 
0 
9 

9999 CALL EXIT 

end 



20021 

6006 

6007 

* 

* 

“V 

* 

* 12 

* 

* 8 

* 9 

* 6 




Cards added to original program 

Modify these cards according to particular geometry considered 


5‘JBROUT INESTrf:- AM ( NPRT , H J', MK , ICYCLE . MCYCLE . HCOVJJ. OX. XCtUW.J2tV2.w2, 
IT) 

DIMENSION Dn.0* < 21, 20 ) , CTDX( 21 , 20 ) , RPSAy (21 ,20 ) , TnS 4V { 21 , 20 ) , 

1RR(?1«20) ,U2( 21,20 ) ,y2( 21,20 ) ,U2< 21,20 ) ,T(2Q) ,C?(?i ,20) ,£'< <21,20 
DIMENSION UJK(21,20) , F JKA ( 21 , 20 ),, FJKN{21,20), QhDOT ( 21 , 20 ) , 

1 QF RAM AX (21,20) , OFRA -*W (21,20) , D MFRAMY ( 21 * 20 ) .CPSAVC20) ,RRSAV(2Q) , 

2 ?3ARi2l,20) ,<J£lA(21,2C) 

C0i v iM0N/STLlN/S,X2,RR,0RQx,DT0X,RPSAy , THSAV, EM, CP, P-3AR 
F(Ei v i,G)=l,/(l.-KG-l 0 ) *'£>-',E.V2. ) 

Ir ( M COW L* . £0 • 1 « AND* F, CYCLE .£0, 1) 60 TO 202 
IF ( ICVCLE .GT. 1) GO TO 55 
C 

C INITIALIZE at x=xi 

c 

CAX =0 o 
CS =0.0 
CN =0.0 
UR I T£ ( 6 , 33U ) X2 

335 FOR, «1AT(///,5X, -INITIAL STREAMLINES POSITIONS At X= *,£15.5/) 

00 34 K=1,MK 
WRITE(6,6001 ) K,T(K) 

RRSAV(K)=RK(1,K) 

C 3 SAV(K)=CP(1,K) 

DO 33 J=2,MJ,2 

OROxt J,K)=V? (U,K)/U2{ J,K) 

OTDX( J,K)=-U2{ J,K)/(RR( J,K)*U2( J*K) ) 

DTOx ( J ,K ) = 57.295779* OTDX(J«K> 

RPSAV( J,K)=RR( J«K) 

THSAV( J,K)=T(<) 

Y1=RR( J,K)*SIN(T(<)/57.295779) 

Z1=RR(U,K)*C0S(T( 0/57,295779) 

WRITE(G.1S2) J, Y1«Z1,RPSa\H J, K> « DRDX ( J , K ) »RR(J«K),T(K)iQTDX(J»K) 

33 CONTINUE 

34 CONTINUE 
RETURN 

C 

C INITIALIZE ODD NUMBERED STREAMLINES AT X= X2 

C 

202 fc/RI TE ( 6« 335 ) X2 
DO 14 K=1,MK 
WRITE(6,6001) K,T(K) 

00 15 J=1,MJ,2 
. DROX( J,K)=V2( J,K)/U2( J,K) 

DT0X(J»K)=-W2<JtK)/<RR{J,K)*U2(JtK)) 

DTDX(J,K)= 57.295779* DTDX(J,K) 

RPSAy( J,K)=RR( J,K) 

THSAV( J,K)=T(K) 

Y1=RR( J,K)*SIN(T( 0/5 7. 295779) 

Z1 = RR<U,K)*C0S(T< 0/5 7.2957 79) 

WRITE(6.1S2) J, Y1»Z1,RPSAV< J,K) ,DRDX( J,K) ,RR{JtK) ,T(K) ,DTQX( J.’O 

13 CONTINUE - 

14 CONTINUE 
30 TO 211 

35 M TE v i=MK-l 

. I F ( v ; C 0 U W .L E. 0) 30 TO 211 
• MSTART=1 

meno=i 


non nnr> ooo ooo 


GO TO 212 

211 '•l S T A R T = 2 
v i E M D = 2 

212 CONTINUE 

00 79 K = 1 ,MTE V 1 
IF(NPRT.EJ. 1) 

lWRITE(b.OOOl) K,T(K) ^ 

6001 FORMATt/t^X.+K— >*«l2,10x, ! i-TH — < : ,El2.5«/,3x, ::: J * , S X » + Y-.lGXtfcZ - , 5 X , 
l*RPREV*,tSX»*0RDx2i:t0X,+R7RAC-,6x«-TTRAC*,5X»¥ERR0Rl*,5X, *=ERR0R2*/} 
00 7 a J=PS l ART ,|MJ*KEN0 

IF ( RPSAV { J • K ) .GE. 1000. .OR. 7H3AV(J«K) .GE, 1000. ) GO TO 76 

lF<0TDX<J«ft) oGE. lOOOo.OR. DRDX<J*K) .GE. 1000„)SO TO 73 

first Guess of streamline position at x+dx 

RPX2 = RPSAV( J.^)+DROX( J.K)=.-OX 
THX2 = THSAV/< J,K) +0TDX( J,K)*OX 
IF(K.EO.I) GO TO 1001 
IF ( J.EG. 1 ) GO TO 350 

TO LOCATE KK AT WHICH THX2 LIES BETWEEN T ( KK-1 ) ANO T(KK) 


kk=o 

51 KK=KK+1 

I r ( KK .GT. MK } GO TO 66 
I F { THX2-T { KK ) ) 53*53,51 

53 IF(KK .EQ. 1) GO TO 1001 

TO LOCATE JJ AT T(K<-1), aNO RPX2 LIES BETWEEN RR < J J-l , kK- l ) AND RR{, 

KK=KK-1 

JJ=0 

57 JJ=JJ+1 

IF ( JJ .GT. MJJ GO TO 66 
I F ( RPX2 -RR(JJ,KK>) 59,59,57 

59 IF ( JJ,EO« 1) GO TO 66 

TO LOCATE II AT T(KK) , ANO RPXg LIES BETWEEN R.R(II-1,KK) AND RR ( 

KK=KK+1 

II = 0 

60 11=11+1 

I F ( I I .GT. MJ) GO TO S6 
IF(RPX2-RR(II,KK) ) 61,61,50 
61 IF l II .EG.. 1) GO TO s 5 

CML1=ITKX2-T(KK-1> >/(T(KK)-T(KK“l) ) 

Cl'AL2=RPX2-RR( J J-l , KK-1) 

CVAL3=RPX2-RR( II-1,KK) 

V/ALl = vU2( JJ, KK-i) -U2( J J-l, KK-1) > / { RR { J J , KK- 1 ) -RR ( J J- 1 , KK- 1 ) ) 

VAL2=( V2( JJ, KK-1) -V2( J J-l, KK-1) ) / { RR { j J , KK- 1 ) -RR ( J J- 1 , KK- 1 ) 5 
VAL3=< W2( JJ, KK-1) -*2< JJ-1, KK-1) )/{RR( JJ,KK-1 ) -RR( J J-l, KK-1) ) 
V\/AL1={U2{1I,KK)-J2(II-1,KK) ) / ( R R { 1 1 , KK ) -RR ( 1 1 -1 , KK } ) 

. V\/AL2= ( V2 ( 1 1 , KK ) - V2 ( 1 1 -1 , KK ) ) / ( RR ( 1 1 , KK ) -RR ( 1 1 -1 , KK } ) 

V V' A L 3 = ( W2 ( 1 1 . KK ) - W2 { I I-i , KK > ) / ( HR (T I , KK ) -RR ( 1 1 -1 , - K } 1 

U2X 2= ( U2( JJ-1,KK-1)+VAL1*C\/AL2) *(1.-CVAL1) + (U2( 11-1, Xis)-- v'V'ALi- 

.1C\/AL3) + CVAL1 

V2X2=(U2( J J-l, KK-1 )+VAL2*CVAL2) s < 1 . -CVAL1 ) + ( V2( II ,-•'»/ AL2' 

1 CV AL3 ) +CVAl_l 


non 


y 2 X 2 =(W 2 ( JJ-l»KK-l)-5-VAL5vCvAL2)«(l.-CvALl) + (W2{II-l»K^H-iy\/AL3* 

X CVAL3)^CUAL1 
350 Ir<J.NL.l) SO TO 360 

U2X?=U2{l.K)-r(U2< itK+l)-U2(i»K) ) / ( T ( K + l ) -T ( K ) ) * C ThX 2- T ( K ) ) 

U2X2 = \/2(1,K) + ( U2( 1«K+1)-V2'(1«K) ) / ( T ( K+l ) -T ( X > )*{T-X2-7 ( <) ) 
l.j2X?=ij2(l,K) + ( W2( 1«K + 1)-W2( 1*K) )/ (T{K + 1 ) - T ( K i )=MTHX2-T ( <) ) 

3&0 0R0X2=V2X2/U2X2 

0TDX2=-W2X2/( ’TPX2SU2X2) 

0 T0X2 = 57,255779* OTOX2 

RTRAC=RPSAV/ ( J* K ) +0 . 5*OX> ( DrDX ( J, K ) +DRQX2 > 

TTRAC=THSAV { J,K) -0 . 5*0X* ( 0T0X { J,K)+DTDX2) 

T£STl=AeS(KTRAC-R J X2 ) /RTRAC 

TEST2=ABS [ T TRAC-’T-SX2 >/TTRAC 

Y1=RTRAC»S1N( TTRAC/57, 295779) 

Zl=RTRAC*COS( ITRAC/57. 235779) 

GO TO 62 

PLANE OF SYMMETRY 

1001 KK = 1 
11 = 0 

IF(J.EQ.I) GO TO 370 
1060 11=11+1 

IF ( I I ,GT , M J ) GO T3 66 
IF(RPX2-RR(II,KK) ) 1061,1061,1060 
1061 IF( II.EQ.DG0 TO 56 

C\/AL3=RPX2-RR( II-1*KK) 

VV/AL1=(U2<1I*<K)-J2( II-1,KK) ) / ( RR ( 1 1 , KK ) -RR { 1 1 -1 , KK ) ) 
VVAL3={W2(II,<K)-*I2(II-1,KK ) ) / ( RR ( II , KK ) -RR ( 1 1 -1 , KK ) ) 

U2X2=U2< II-l,<K)+i/v/ALl*CVAL3 
V2X2= V2(II-1»KK) + VVAL2«CVAL3 
DRDX2=V2X2/U2X2 

370 IF(J,£q,1) 0R’JX2= </2 ( 1 , K ) /U2 ( 1 , K ) - — 

0T0X2 = 0.0 

RTRAC=RPSAV( J,K) +0.5*3X*(DR0X{ J,K)+0RDX2) 
ttrac=o.o 

T-ST1=ABS(KTRAC-R 3 X2 )/RTRAC 

TEST2 = 0.0 
Y 1 = 0 . 

21 = RTRAC 
62 CONTINUE 

GO TO ( 60 0 6 , 6007) NPRT 

6006 tdRI TE ( & , 162 } J,Y1,Z1, RPSAV ( J , K ), DRDX2 , RTRAC , TTRAC , TEST1 , TEST2 
162 F0R v .AT(5X,15,9El2.5) 

6007 CONTINUE 
GO TO 68 

66 WRITE<6,67) K< « I I ♦ J J » J » K , RPX2 , IHX2 » RPSAV ( J , K ) , X2 

67 FORMAT!* CAM NOT LOCATE STREAM LINE AT THIS PT . * , 51 = ♦ £15 . 5 ) 
RTRAC=10GQ. 

ttrac=igoo. • 

0R0X2=10G0. 

OTDX2 = 1000. . _ - - 

68 CONTINUE 

RPSAV (J,K)=RTRAC 
THSAV( J,K)=TT+AC 
. 0R0X( J,K)=0R0X2- 

OTOX( J,K)=0T0X2 
78 CONTINUE 


non 


C 

c 

c 


79 CON TIMUE 

calculate incremental BODY FORCE AND MOMENT 
Ir(MPRT.EQol) UR I rE(S, 312) X2 

312 F0R'-'.AT(///.5X,*S0JY FORCE COMPONENTS AND MOMENTS UP TO X = ;- , E15 . 5 , / 
1 1 0 X i *K *«5X,ftC AxIAL*, 5x» a-CNORtTAL* , 5X * *CSI0£*5X« *C::Y=5= « / ) 

mjmi=mj-i 
M K M 1 = M K - 1 
00 315 K = 1 » MK “11 

R1AVG={RRSA\/{X->1) 4-RRSAV(K) ) /2 e 
R2AVG=(RR(1,K*1)-HR(1,K) )/ 2 o 
RAVG=(R1A^/i3vR2Av/S)/2. 

DRAv/G=R2A\/G-r1AVG 
OT = T ( K+l ) -T ( K ) 

DSAX = RAV/G-icDF^DRAVG /57.29S779 
f {QYA\/G=(DY2 + 0Yl)/2. 

| i)Y2sRR(l,<+l)*Sl,M( T(K + l)/57o 2957 79 ) -RR( 1 , K ) *SIN ( T ( K) /57. 295779) 
LOYl = RRSAV(K + l) v'SI\!<T(’< + 1 5/57.295779) -RRSAV(K)-SlN(T(K)/57.295779) 
QSN r DYAl'b-tDX 

DZ2=RR(l t K+l)*COS( T(K+1) /57 o 295779)-RR(i»K)*C0S(T{K)/57. 295779) 

OZ1 = RRSaV«K- 1 '1) ^'C03( T( K + l) /S7.295 779 )-RRSAV(K) *COS ( T { x ) /37 . 295779 ) 
0ZA\/G=(D22 + DZl)/2. 

0SS=D2AVGiDX 

CPAVG = {CPSAy(X)+C p SA\/(K + l)+CP(l,i<)+CP(l < K+l) l/tf. 

ocax=cpa\/g*osax 

OCN=£PAVG*USN 

OCS=CPAvG^uSS 

RRSAV (K)=RR(1,K) 

CPSAV(K)=CP(1,K) 

T3AR=(T(K+1)+T(K) )/2. 

R3ARSR1AVG+0RAV&/2, 
y3AR=RBAR*SlN(TBAR/57 
Z3AR=RBAR^C0S( T3AR 
X3AR=X2+OX/2. 

DCMY=XBAR*OCN -ZBAR^DCAX 
C AX =CAX +OCAK 

CM =CN + 0CN 

CS =CS +0CS 

CrtY =CMY +QCMY 

IF(NPRT.EQ.I) W R I f E ( 6 1 3 2 0 ) 

320 F0RMAT(5X«15,‘fE15oH) 

313 CONTINUE' 


295779) 

/57. 295779) 


K,CAX«CN, CStCMY 


calculate incremental mass flow Rate, Ram force, pitching moment 


IF ( NPRT . NE« 1 ) RETURN 
00 301 K=1«MK 
00 302 J=1,MJ 

S0=SQRT ll.+(V2( J,X)**2+W2( J, K ) =s*2 ) /U2 { J«K)**2) 

0JK( J,K)=EM{ J,K)/SGRT ( F ( EM ( J , K ) , S ) ) *P3AR ( J , K ) /SO 
FJ*A(J»K)= ( l.+G*£M( J c K ) ^*2/SQ/S2 ) *?3AR ( J,K)-1« -- - 

FJKN(J,K)= (G^£M( J,K)^'v 2)/S3/5Q?P3AR{ J,K)*<V2( J , K ) -COS ( T ( <) /57 « 25: 
1779 ) -w2 ( J , K ) *SlN { T ( K ) /57 . 295779 ) )/U2( J,K) 

302 CONTINUE 
301 CONTINUE 

00 303 K=1 ♦ MKM1 
SJM1=0.0 


S~s 

SUM 2=0.0 
SUH3=0»Q 
SJW*4 = 0.0 
S U M 5 = 0 . 0 
DO 304 J=l»MJ l U 

OJKaVG = (OJfUJ.K) -i- 0 JK ( J + l » K ) + 0 JK ( J ♦ K-5-1 > + Q JK { J+i , (<-*-1 ) ) / '4 . 
FJXflAV=( FJ*A ( J'K) <-FJKA( J-M , r< ) +F J<A < J . K+l ) -*-F JKA ( J'-l t <-5-1 ) )/!», 

F J K M A V = < F JKN < J ♦ K ) f F JKN < J + l , K j +F JXN ( J , K + l ) +F JKM ( Jvl , K+i } ) /'+ . 

R1AVG = {RK( J»K)+RR( J,K+1) >/2„ 

R2A\/G=(RR ( J + l «K)+RR< J + l, K + l > >/2. 

RA\/5=(RlAVbvR2AVS)/2. 

DRAi /5= R2AVS-R1A 1/3 
DT=T(K*X)-1 (K) 

OAJKsRAVG^OT^'JRAVo/ 57, 295779 
SlMl=SUdl+UAd« 

0ELA ( J«K) =bUMl 
SJM2=SU*-i2+UJKfM/G* JAJK 
D u iDOT(0,K)=sU v ''2 
OrRAMA=FJKAAV*OAJ< 

SUM3=SUM3+0FRAtfA 
OFRftKAXl J»rt)=SUfl3 
QFR A MNO=F Jrt M A V/ * D A JK 
SJM4=SUM4->-0 FrAM.\I 0 
OFRAM'J( J,K)=SUM4 
R3AR=RiAV3+0RAVG/2. 

T3AR=(T(K + l)-t-r (K) )/2. 

Y3AR = R3AK*SIN( T 3AR/57 „ 295779 ) 

Z3AR=RBAR*C0S( T3AR/57. 295779) 

SUM5=SUM5+X2* OF RAMt\10-ZBAR*0 rR AM A 
OMFRAMYt U,K)=SUM5 

304 continue: 

303 CONTINUE 1 * 

wRITE(b,300) (K«Ksi«n<Ml) 

300 FORMAT (///*5X» *TA3LE OF INCREMENTAL MASS FLOW RATE* «/. 5X « 12110 ,/ ) 
00 314 J=1»MJM1 

314 WRIT£(6«310) J« (0M007( J,K> «K=1«KKM1> 

310 -0R v tAT(5Xfl5 f 12El0 o 3) 

WRITE(6,330) (K,K=l,fKMl) 

330 FORMAT W//»5X.*INCREMEN7AL AREAS. IN MESH GRID* , / , 5X t 1211 0 t / ) 

00 331 J=1«MJM1 

331 WRITE(6,310) J» (0£LA( JtK) «K=I,rlKMl) 

WRITE(6,303) (K v K=l«K<ni} 

305 F0RMA7(///,5X,*-AXIAL COMPONENT OF RAM FORCE* , / « 5X , 1 21 10 , / ) 

. 00 506 J=1«MJM1 

306 WRITE 1 6 « 310 ) J . ( D = RAMAX ( J , K) ,K=1 « MKhl ) 

WRITE(6,307) ( K » K=1 » f-HMl ) 

307 FOR '•'1AT(///«5X» S - NORMAL COMPONENT OF RAM FORCE* ♦/. 5X , 1211 0 ,/ ) 

00 303 J=1 y M J Ml 

308 WRITE(6,310) J « (OFRAMN ( J, K j , K=1 , MKHl ) 

WRITE(6,509> CK t K = l,i-HMl) 

309 FORMAT ( /// 1 5X » '-MOMENT PRODUCED BT RAM FORCE* ,/. 5X ♦ 12110 ./ ) 

DO 31 1 J=1«MJM1 _ _ 

'311 JRITE(6,310) J » { OMFrAMy ( J « K ) « K=1 « MKMl ) 

BO RETURN 
END 


1 


T A e?L £ jL 




C(voVvx^£c\ I 

Hi W 


INCREMENTAL areas 

1 

1 6.504E 

2 1 . 3 2 4 £ 

3 2.019E 

4 2.73QE 

5 3.4 79 E 

6 4.243E 


7 5.030C 
0 5.B40E 
9 6.672C 

10 7.527E 

11 0.405E 

12 9.305E 

13 1.023E 

14 1.U7E 

15 3 ,23 4F 


~ Q 3 6 
02 1 
-02 2 
-02 2 
02 3 
02 4 
-02 4 
-0 2 5 
-02 6 
-0 2 7 
- 0 ? ft 
-02 9 
-01 1 
-01 1 
-01 3 
1 


1C 1.313E 
17 1.43 5E 
10 1.519C 

19 1.625E-''* f — r -- 

20 3.733E-01 3,/20£-01|3 l 


•01 
• 01 
>01 
>01 


IN MESH 
2 

. 4 5 4 fc- - 0 4 
.315L-02 
, 0 0 4 E « 0 2 
.717002 
.455002 
.23 IE-02 
,992002 
. 794L..02 

.621 L-0? 
.470 L -02 
o 3 4 3. E - 0 2 
, 2 3 5 1- - 0 2 
.015001 
. 109001 
,205001 
. 304L-0 l 
,404001 
.^0 7001 
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1 .925C-02 
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1 .911E-05 
2.59200 2 
3 . 294002 
4 o 010E-0 2 
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U 
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13 
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, 7500 i 00 
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,‘4 860:00 
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1,0490+00 
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1. 7260+00 
1.9590+00 
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47 U +00 


0. 
1. 
1. 
1 0 
J . 
1. 
2 o 
?. 


? . 72.20 + 00 
2.976+ ; 00 
3,2550)00 
5.497L+00 
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9.957E-0 L 

10 

1.32UC+00 

J , 2 3 2 1. } 0 0 

1.04500 0 

19 

1 , 380E f-00 
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J , 34U * oif| l , J 47E*0 0 


4,7240 
9 . 462 E- 
i,4?oe.. 
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, 27.1C 


-02 9 
0 ? 1 
02 2 


-01 
-01 
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